Low temperature properties of the infinite-dimensional attractive Hubbard model 



Department of Physics, 



Akihisa Koga 
Tokyo Institute of Technology, 



Tokyo 152-8551, Japan 



Philipp Werner 
Theoretische Physik, ETH Zurich, 8093 Zurich, Switzerland 
(Dated: June 24, 2011) 

We investigate the attractive Hubbard model in infinite spatial dimensions by combining dynam- 
ical mean-field theory with a strong-coupling continuous-time quantum Monte Carlo method. By 
calculating the superfluid order parameter and the density of states, we discuss the stability of the 
superfluid state. In the intermediate coupling region above the critical temperature, the density of 
states exhibits a heavy fermion behavior with a quasi-particle peak in the dense system, while a dip 
structure appears in the dilute system. The formation of the superfluid gap is also addressed. 



I. INTRODUCTION 

Ultracold atomic systems have attracted a lot of in- 
terest following the demonstration of Bose-Einstein con- 
densation (BEC) in a Rb atom system Among the 
many interesting issues, the nature of the superfluid state 
in fermionic systems [2| is a widely studied topic, since 
a crossover between a weak-coupling BCS-type super- 
fluid state and a strong-coupling BEC-type superfluid 
state (BCS-BEC crossover) has been observed in experi- 
ments Recently, the appearance of a pseudogap in 
the density of states has also been suggested @, which 
stimulates further theoretical investigations on the super- 
fluid state and related phenomena. 

Fermi gas systems have been studied theoretically in 
much detail and it has been clarified that a pseudogap 
phenomenon indeed appears in the BCS-BEC crossover 
region above the critical temperature More re- 

cently, the gap structure in the superfluid state has been 
discussed [TJ [T^]. It has been pointed out that upon 
decreasing temperature, the pseudogap disappears at a 
certain temperature inside the superfluid phase and the 
superfluid gap opens instead fl2l ]. On the other hand, 
in fermionic optical lattice systems described by the at- 
tractive Hubbard model, such dynamical properties at 
finite temperatures have not been discussed. It is also 
important to clarify how the particle density as well as 
the interaction strength affect the low energy properties 
in the pseudogap region. Therefore, it is useful to study 
the stability of the superfluid state systematically and to 
clarify how the gap structure appears in the density of 
states at low temperatures. 

In this paper, we consider the infinite-dimensional 
attractive Hubbard model to discuss its low temper- 
ature properties. To take correlation effects precisely 
into account, we combine dynamical mean-field theory 
(DMFT) (ill - fl(| with a continuous-time quantum Monte 
Carlo (CTQMC) method [13 • This unbiased technique 
enables us to discuss the stability of the superfluid state 
quantitatively. By tuning the particle density, interaction 
strength, and temperature, we determine phase diagrams 
of the system and discuss the formation of the gap in the 



density of states. 

The paper is organized as follows. In §2, we introduce 
the attractive Hubbard model and briefly summarize our 
theoretical approach. We discuss how the superfluid state 
is realized at low temperatures in §3. In §4, we focus on 
the dilute system in the BCS-BEC crossover region. We 
clarify how the gap structure appears in the density of 
states and discuss the difference of low energy properties 
in lattice and Fermi gas systems. A brief summary is 
given in the last section. 



II. MODEL AND METHOD 

We consider two-component fermions in an optical lat- 
tice, which may be described by the following Hubbard 
Hamiltonian, 



n = -t 4* c ]<? - U Y1 n i j t n ih 

(ij),<r i 



(1) 



where Ci a (cj^) is an annihilation (creation) operator 
of a fermion on the ith site with spin cr(=t,4-), and 
rii a = c\ a Ci a . U is the onsite attractive interaction and 
t is the transfer integral between sites. The low-energy 
properties of the attractive Hubbard model have been 
studied in one dimension 16 - 23||. two dimensions [24l - 
|2(1 | and infinite dimensions 27 435| . It is known that 
the superfluid ground state is always realized in two and 
higher dimensions, where the BCS-BEC crossover has 
been studied in detail. However, dynamical properties 
have not been studied yet in detail in the intermediate 
correlation and temperature region. In particular, the 
question whether and how a pseudogap appears in the 
density of states above the critical temperature, has not 
been answered. In this paper, we systematically inves- 
tigate the low temperature properties in the attractive 
Hubbard model by varying the particle density, interac- 
tion strength, and temperature. We then clarify how the 
gap structure appears in the density of states. 

For this purpose, we make use of DMFT [ilj - flfjl ]. In 
DMFT, the original lattice model is mapped to an ef- 
fective impurity model, where local particle correlations 
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are accurately taken into account. The lattice Green's 
function is obtained via a self-consistency condition im- 
posed on the impurity problem. This treatment is exact 
in infinite dimensions, and the DMFT method has suc- 
cessfully been applied to strongly correlated fcrmion sys- 
tems. In DMFT, we take into account dynamical correla- 
tions through the frequency-dependent self-energy. This 
allows us to discuss the stability of the s-wave superfluid 
state more quantitatively than in the static BCS mean- 
field theory. 

The lattice Green's function is given by 

G -1 (fc, iu> n ) = iu> n fr Q + (/i - e fe ) <7 Z - £ (ioj n ) , (2) 

where fi is the chemical potential, &o and a z are the iden- 
tity matrix and the z-componcnt of the Pauli matrix, 
is the dispersion relation for the non-interacting system, 
and ui n is the Matsubara frequency. S(«w n ) is the local 
self-energy in the Nambu formalism. The local lattice 
Green's function is obtained as 



(3) 



G, which are given 32[ by 



G(iuj n ) — / dkG(k,iu; n ). 



In this paper, we use a semi-circular density of states, 
p{x) = 2/nDy/l - (x/D) 2 , where D is the half band- 
width, which corresponds to an infinite coordination 
Bethe lattice. The self-consistency equation (3(| is then 
given by 
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There are various numerical methods to solve the 
effective impurity problem. To study the attractive 
Hubbard model systematically, an unbiased and ac- 
curate numerical solver is necessary, such as the ex- 
act diagonalizatio n[32l 34 . l37l or the numerical renor- 
malization group particularly powerful 

method for exploring finite temperature properties is 
CTQMC [l7|. In our previous paper 35|, we have used 
the CTQMC method in the continuous-time auxiliary 
field formulation [4l[ to study the imbalanced attrac- 
tive Hubbard model. However, this is a weak-coupling 
approach, which is not suitable for a systematic inves- 
tigation of the low-temperature properties in the strong 
coupling regime. Hence, our previous discussion was re- 
stricted to the weak coupling region. Here, we employ the 
complementary strong coupling version of the CTQMC 
method (42[, which is more efficient in the large \U\ re- 
gion. We use this method to investigate the attractive 
Hubbard model both in the weak and strong coupling 
regimes. Some details of the implementation are ex- 
plained in the appendix. 

In this paper, we use the half bandwidth D as the 
unit of energy. We then calculate the pair potential A, 
double occupancy d, internal energy E, and specific heat 



A = (c n c ti ) = G 12 (r = 0+), 
d = (n lt n 4 }, 



(5) 
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In addition to these static quantities, we deduce the 
density of states by applying the maximum entropy 
method [43T - |45| to the Green's function. We then dis- 
cuss how the gap structure appears in the system. 



III. STABILITY OF THE SUPERFLUID STATE 

We first consider the attractive Hubbard model with 
different band fillings to discuss how the superfluid state 
is realized at low temperatures [3ll - l34^ . Figure [1] shows 
the pair potential A and the double occupancy d at a 
fixed temperature T = 0.02. In the noninteracting case 
(U = 0), a normal metallic state is realized, with A = 
and d = n 2 . Attractive interactions lead to the formation 
of Cooper pairs and an increase in the double occupancy, 
as shown in Fig. [T](b). More precisely, at a certain criti- 
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FIG. 1: Pair potential (a) and double occupancy (b) as a 
function of the interaction at temperature T = 0.02. Circles, 
squares, and triangles show the results for the system with 
n = 0.5, 0.25 and 0.125, respectively. 
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cal interaction U c i, a phase transition occurs to a super- 
fluid state and the pair potential is induced, as shown in 
Fig. [T] (a). A cusp singularity appears in the curve of the 
double occupancy although it is not visible on this scale. 

To clarify how the phase transition affects dynamical 
properties, we show the density of states deduced by the 
MEM in Fig. H When U < U cl , the system is metallic 
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Eventually, the pair potential vanishes and a phase tran- 
sition occurs to the normal state at another critical point 
U C 2, as shown in Fig.[T](a). In contrast to the BCS region, 
a large gap remains in the density of states and the phase 
transition little affects its features, as shown in Fig. [2] 
This originates from the fact that in the region (U > 5), 
paired bosons exist in the normal state as well as in the 
supcrfluid state. Therefore, we conclude that a BEC-type 
supcrfluid state, which can be regarded as the condensa- 
tion of paired bosons, is stabilized at low temperatures. 
The critical values are obtained as U C 2 ~ 12.5, 11.4 and 
9.7 for n = 0.5, 0.25 and 0.125. 

In the superfhhd phase, the BCS-type state (U ~ 1) is 
adiabatically connected to the BEC-type one (U > 5). A 
BCS-BEC crossover thus occurs between the two states. 
Figure [3] shows the temperature dependence of the pair 
potential for three cases with U = 1,2 and 5. When 
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FIG. 3: Pair potential as a function of temperature in systems 
with U = 1.0 (a), 2.0 (b), and 5.0 (c). Circles, squares, and 
triangles represent the results for the systems with n = 0.5, 
0.25 and 0.125, respectively. 
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FIG. 2: Density of states for the systems with n = 0.5 (a), 
0.25 (b) and 0.125 (c) at T = 0.02. 



with a finite density of states around the Fermi level. If 
the system enters the supcrfluid state, weakly coupled 
Cooper pairs are formed, which give rise to a tiny su- 
pcrfluid gap together with coherence peaks at its edges. 
This characteristic behavior is found for all the band fill- 
ings considered, so we can say that a BCS-type supcr- 
fluid state is realized in the region (U ~ 1). By examin- 
ing the critical behavior of the pair potential, we deduce 
lid ~ 0.62,0.64 and 0.68 for n = 0.5,0.25 and 0.125, 
respectively. 

When U ~ 5, the pair potential has a maximum and 
the double occupancy approaches the particle density 
[d ~ n), as shown in Fig. [TJ This means that most 
fermions are tightly coupled to form paired bosons, which 
stabilizes the superfluid ground state. In this case, the 
energy scale should be t 2 /U, characteristic of the hopping 
for the paired bosons. Therefore, a further increase in the 
attractive interaction effectively increases the tempera- 
ture of the system, making the supcrfluid state unstable. 



the temperature is decreased, the pair potential is in- 
duced at the critical temperature, where the phase transi- 
tion occurs from the normal state to the supcrfluid state. 
We find that the pair potential increases monotonically 
with decreasing temperature and saturates at low tem- 
peratures (T < 0.01). Extrapolating the pair potential 
to zero temperature At=o, we also deduce the quantity 
q = UAt=o/T c , as shown in Fig. [4] In the weak cou- 
pling limit (U — > 0), a BCS-type superfluid state is re- 
alized, where this quantity is independent of the band 
filling and takes the universal value qo = 1.764, accord- 
ing to the simple BCS mean-field theory. Increasing the 
interaction strength, q monotonically increases, in a way 
which depends on the band filling. In the strong cou- 
pling region, the quantity is proportional to the square 
of the interaction strength q ~ U since the pair poten- 
tial At=o should saturate and T c ~ 1/U. We note that 
these results differ considerably from those obtained by 
the simple BCS mean-field theory, where quantum fluc- 
tuations are not taken into account properly. In fact, the 
critical temperature is always overestimated (as shown 
in Fig. [S]) and qscs is almost constant, as shown by the 
dashed lines in Fig. 2] Therefore, it is crucial to take into 
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critical temperature in systems with U = 1,2 and 5. 
It is found that the critical temperature is slightly de- 
creased away from half filling. When the particle density 
approaches the dilute limit (n — > 0), the critical tempera- 
ture rapidly decreases, as T c ~ n a , where a = 0.2 ~ 0.3. 

In this section, we have discussed how the superfluid 
state is realized at low temperatures and determined the 
phase diagram of the attractive Hubbard model. In the 
following section, we focus on the dilute system to ex- 
amine dynamical properties in the BCS-BEC crossover 
region. We then clarify how the gap structure appears in 
the density of states. 



FIG. 4: The ratio q/qo as a function of the interaction U. 
qo is the universal value in the weak-coupling limit obtained 
from the BCS theory. Circles, squares, and triangles repre- 
sent the results for the system with n = 0.5, 0.25 and 0.125, 
respectively. 



account dynamical correlations correctly in the attractive 
Hubbard model with U > 1. 

By performing similar calculations, we determined the 
phase diagram shown in Fig. [SJ In the weak coupling 
region, the BCS-type superfluid state is realized. On the 
other hand, the BEC-type superfluid state is realized in 
the strong coupling region, where the phase boundary 
scales with the inverse of the interaction, as discussed 
above. The critical temperature reaches its maximum 
value in all cases n = 0.5,0.25 and 0.125 in the inter- 
mediate region (U ~ 2), where the BCS-BEC crossover 
occurs. We also find that a decrease of the particle den- 
sity monotonically shifts the phase boundary. The inset 
of Fig. [5] shows the particle density dependence of the 
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FIG. 5: Phase diagrams for the attractive Hubbard model. 
Circles, squares, and triangles indicate the phase boundaries 
in the systems with n — 0.5,0.25 and 0.125. Solid (open) 
symbols are the results obtained by DMFT with CTQMC 
(the BCS theory). Lines are guides to eyes. Inset shows that 
the critical temperature as a function of the particle density. 
Open circles, squares, and triangles are the results for [7 = 1,2 
and 5. 



IV. LOW-TEMPERATURE PROPERTIES NEAR 
THE BCS-BEC CROSSOVER 

In the section, we study low temperature properties in 
the intermediate coupling region, where the BCS-BEC 
crossover occurs. It is known that in the dense system, 
the BCS-BEC crossover occurs at an interaction which 
is slightly smaller than the value corresponding to the 
maximum and no pscudogap appears near the 

phase boundary 35|. On the other hand, it has been 
reported that a pseudogap always appears at the critical 
temperature in the Fermi gas system (T^]. Therefore, 
it is necessary to clarify how the pseudogap behavior is 
realized in the dilute system. To clarify this, we focus on 
the attractive Hubbard model with a low particle density 
(n = 0.125) and compute spectral functions in the BCS- 
BEC crossover region. 

We note that the gap structure appears even in the 
normal state if the preformed pairs are stabilized by the 
large attractive interaction. This means that it is not 
directly related to the realization of the superfluid state, 
but rather to the crossover between the metallic and in- 
sulating state, which occurs in the normal phase. To 
make this clear, we first examine the interaction depen- 
dence of dynamical properties at the high temperature 
T = 0.087. Figure [6] shows the density of states for the 
system with n = 0.125. We clearly find double peaks 
in the case U = 1.5, where the metallic state is real- 
ized. A remarkable point is that a large density of states 
appears around the Fermi level. This is reminiscent of 
heavy fermion behavior in the repulsive Hubbard model, 
where a quasi-particle peak develops in the metallic state 
close to the Mott transition. In the attractive case, if the 
ground state is restricted to be paramagnetic, a pair- 
ing transition occurs between the metallic and insulating 
states at any filling f30l ] . The results in the paramagnetic 
state at a lower temperature (T = 0.02) are shown as 
the dashed lines in Fig. [B] It is found that the sharp 
quasi-particle peak grows in the region U < 2. There- 
fore, we conclude that the large density of states near 
the Fermi level, which appears above the critical tem- 
perature at U — 1.5, results from particle correlations. 
We wish to note that this behavior is characteristic of 
the lattice model, in contrast to the interacting Fermi 
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FIG. 6: Solid lines indicate the density of states at T = 0.087 
in the system with n — 0.125. Dashed lines indicate the 
results at T — 0.02 in a system which is restricted to be 
paramagnetic. 
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FIG. 7: The pair potential (a), the double occupancy (b), the 
internal energy (c), and the specific heat (d) as a function of 
temperature around the BCS-BEC crossover. 



gas system where a gap behavior appears at the critical 
temperature [l2j . 

Further increase in the attractive interaction leads to 
the crossover from the heavy metallic state to the in- 
sulating state. It is found that the quasi-particle peak 
collapses and a dip structure appears instead around the 
maximum T c value (U ~ 2.0). Its width continuously 
grows with increasing interaction, as shown in Fig. [fj] 
Thus, strong pairing correlations stabilize the gap struc- 
ture in the vicinity of the Fermi level. 

In the intermediate coupling region the crossover be- 
tween the metallic and insulating states occurs, which is 
associated with the pairing transition in the normal state, 
as discussed above. It is known that at zero temperature, 
this transition point is shifted away from half filling, e.g. 
U c ~ 3.0^.4 and 1.12 in the cases with n = 0.5,0.25 
and 0+ [301 . Therefore, we can say that the low energy 
properties around the critical temperature depend on the 
particle density as well as the interaction strength. This 
should have important implications for the observation 
of the pseudogap behavior in fermionic optical lattices, 
where these parameters can be controlled experimentally. 
Namely, in the BCS-BEC crossover region, heavy fermion 
behavior appears in the dense system |35| . while a dip 
structure (pseudogap behavior) is found in the dilute sys- 
tem. 

Let us consider the temperature dependence in the 
dilute system in the BCS-BEC crossover region. We 
calculated the pair potential, double occupancy, inter- 
nal energy, and the specific heat, as shown in Fig. [7] 
At the critical temperature, the pair potential is in- 
duced, and a cusp singularity appears in the curves of 
the double occupancy and the internal energy. A jump 
singularity appears in the curve of the specific heat, 
and its height increases with increase of the interaction. 
We note that the cusp singularity in the double occu- 
pancy vanishes at a certain interaction ?7*(~ 1.75), i.e. 



Ad/AT\ T 

=t c +s — Ad/ AT\t=t c -8- This is due to the 
crossover between the heavy metallic state and the in- 
sulating state. When the temperature is decreased in 
the case U < U*, the double occupancy reaches a maxi- 
mum at a certain temperature T max and is decreased as 
temperature is lowered to the critical temperature. This 
suggests the formation of the heavy fermion state in the 
region with T c < T < T max . In fact, we find an en- 
hancement of the quasi-particle peak in Fig. [5] (a) . On 
the other hand, when U > U*, the attractive interac- 
tion stabilizes preformed pairs at high temperatures and 
thereby the decrease in the temperature simply increases 
the double occupancy. Thus the dip structure in the den- 
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FIG. 8: Density of states in the dilute system (n = 0.125) 
with U = 1.5 (a) and 2.0 (b), where T c ~ 0.076 and 0.083, 
respectively. 
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sity of states is stabilized above the critical temperature, 
as shown in Fig. 0(b). 

Below the critical temperature, the static quantities 
change gradually, as shown in Fig. [7] When U = 1.5, the 
decrease of the temperature increases the pair potential, 
where the heavy quasi-particle peak collapses and the 
supcrfluid gap opens around the Fermi level, as shown 
in Fig. [8] (a). On the other hand, when U = 2.0, the 
dip structure found at the critical temperature smoothly 
evolves into the gap structure at low temperatures. This 
behavior is also observed in the system with a lower par- 
ticle density n = 0.025 at the temperature T/T c > 0.17. 
Therefore, we can say that the region with a dip struc- 
ture at high temperatures is adiabatically connected to 
that with the supcrfluid gap at low temperatures. This 
is in contrast to the results for the three-dimensional 
Fermi gas system [l2|]. In the gas system, decreasing 
the temperature around the BCS-BEC crossover smears 
out the pseudogap structure, and the superfluid gap de- 
velops below a certain temperature. This may suggest 
that fc-dependent correlations plays a crucial role in re- 
alizing these dynamical properties. This topic is beyond 
the scope of our paper, but it is important to clarify how 
two kinds of gap structures appear in a low dimensional 
optical lattice system, by taking into account spatial cor- 
relations as well as dynamical correlations. 



SUMMARY 



We have investigated the attractive Hubbard model 
in infinite spatial dimensions by combining dynamical 
mean-field theory with a strong-coupling continuous-time 
quantum Monte Carlo method. By calculating the su- 
pcrfluid order parameter and the double occupancy, we 
have systematically studied the stability of the super- 
fluid state and determined the phase diagram of the sys- 
tem. By computing the density of states, we have found 
that the gap structure is strongly affected by the inter- 
action strength and the particle density, which is asso- 
ciated with the pairing transition in the normal state. 
Namely, around the BCS-BEC crossover, a dip structure 
appears in the dilute system while heavy fermion behav- 
ior is found in the dense system. We have also examined 
the dynamical properties in the superfluid state and have 
clarified that the dip structure above the critical temper- 
ature continuously evolves into the supcrfluid gap with 
decreasing temperature. 

In this paper, we have considered the dynamical prop- 
erties characteristic of the infinite-dimensional lattice 
model, which are qualitatively different from those in the 
Fermi gas system. It is an interesting problem to clarify 
how the low energy properties are changed by adding the 
lattice potential to the Fermi gas system, which is now 
under consideration. 
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Appendix 

In this study, we have used the strong-coupling version 
of the CTQMC method to solve the effective Anderson 
impurity model [42| ■ This method is based on a stochas- 
tic sampling of an expansion of the partition function 
in powers of the impurity-bath hybridization. A given 
Monte Carlo configuration can be represented by a col- 
lection of segments {r^, r icr }, where r z ' CT (ri CT ) is the start- 
ing (end) time for the i-th segment with spin a. These 
segments represent time intervals in which an electron of 
spin a resides on the impurity. The Monte Carlo simula- 
tion proceeds via local updates, such as insertion/removal 
of a segment or empty space between segments (called 
anti-segment), or shifts of segment end-points. 

Here, we discuss updates which improve the sampling 
efficiency in the strong-coupling region. For \U\ > 3, 
the large energy cost for inserting or removing (anti-) 
segments leads to a high rejection rate for proposed in- 
sertion/removal updates. In fact, the local (impurity) 
contribution to the Monte Carlo weight is given by 



Wlc 



cxp 



-E 



Ul 



overlap 



(9) 



where Ef is the energy level at the impurity site and 
l a and loveriap are the total lengths of a segments, and 
the overlap between up and down segments. Thus, in 
the strong-coupling regime, at low temperature, the ac- 
ceptance probability for a generic update will be expo- 
nentially suppressed. One possibility is to take this ex- 
ponential dependence into account on the level of the 
proposal probabilities. Another possibility to overcome 
this bottleneck is to consider updates which change the 
configuration for both spins simultaneously. When both 
spins arc flipped between occupied and unoccupied states 
in a certain time interval, the energy change on the impu- 
rity site is not so large. In particular, the corresponding 
energy change is zero at half filling (Ef = U/2). 

We discuss here explicitly one of the simplest of these 
updates, which is schmatically shown in Fig. |H1 We first 
choose a random time interval, defined by a pair of neigh- 
boring creation/annihilation operators {Tj CT , r' i(J \. Note 
that the spin and type (creator/annihilator) is arbitrary, 
so the time interval need not correspond to any segment. 
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The proposal probability is given by p prop = (2iV) _1 
where N is the total number of segments. The update 
consists of flipping both spins in the time interval. Note 
that this process is composed of two standard updates, 
i.e. two insertions/removals of a segment/antisegment, 
or two shift updates. Therefore, it is easy to add these 
updates to an existing CTQMC code. The above updates 
were found to improve the efficiency of Monte Carlo sim- 
ulations in the strong coupling regime. 



FIG. 9: Efficient updates in the strong-coupling regime. The 
two imaginary time intervals correspond to spin up and down. 
Bold lines represent segments, i.e. a time interval, in which 
an electron resides on the impurity, (a) The shaded regions 
indicate examples of time intervals for which an improved 
updated will be proposed, (b) Resulting configuration if the 
update for time intervals 1 and 3 is accepted. These updates 
correspond to two shifts of segment end-points, (c) Resulting 
configuration if the update for intervals 2 and 4 is accepted. 
Here, the update corresponds to two insertions/removals of a 
segment /antisegment. 
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